FTICR Visualizations

# Clean workspace
rm(list=ls())

library("ftmsRanalysis")
library("vegan")
library("ggplot2")
library("ape")
library("ggVennDiagram")
library("reshape2")

Load Peak Object

load("ficusPeakObj.Robj")
str(peakObj)

List of 3
 $ e_data:'data.frame': 3704 obs. of  33 variables:
  ..$ Mass: chr [1:3704] "200.0564425" "200.092833" "201.004077" "201.0404701" ...
  ..$ H1  : num [1:3704] 21.7 20.2 20.2 22.3 NA ...
  ..$ H2  : num [1:3704] NA NA NA 20.3 NA ...
  ..$ H3  : num [1:3704] 21.3 NA NA 22.2 NA ...
  ..$ H4  : num [1:3704] NA NA NA 20.9 NA ...
  ..$ H5  : num [1:3704] 21.1 NA NA 21.8 20.3 ...
  ..$ H6  : num [1:3704] 20.3 NA 20.2 21.2 NA ...
  ..$ H7  : num [1:3704] 21.1 NA NA 21.5 NA ...
  ..$ H8  : num [1:3704] NA NA NA 20.7 NA ...
  ..$ H9  : num [1:3704] 21.2 20.6 20.3 22.7 NA ...
  ..$ H10 : num [1:3704] NA NA NA NA NA ...
  ..$ H11 : num [1:3704] 21.7 20.7 NA 22.2 20.6 ...
  ..$ H12 : num [1:3704] 20.4 NA NA 21.2 NA ...
  ..$ H13 : num [1:3704] 21.5 NA NA 21.8 NA ...
  ..$ H14 : num [1:3704] 20.8 NA 20.2 21.9 NA ...
  ..$ H15 : num [1:3704] 21.8 NA 20.7 22.6 21 ...
  ..$ H16 : num [1:3704] 21 NA NA 21.8 NA ...
  ..$ C1  : num [1:3704] 21.7 20.2 20.4 22.8 NA ...
  ..$ C2  : num [1:3704] NA NA NA 21 NA ...
  ..$ C3  : num [1:3704] 21.1 20.3 20.2 22.4 20.8 ...
  ..$ C4  : num [1:3704] 20.2 NA NA 21.7 NA ...
  ..$ C5  : num [1:3704] 21.2 NA 20.2 22.7 20.3 ...
  ..$ C6  : num [1:3704] 20.3 NA 20.4 21.7 NA ...
  ..$ C7  : num [1:3704] 21.9 NA 20.8 23 20.3 ...
  ..$ C8  : num [1:3704] 20.3 20.2 NA 21.5 NA ...
  ..$ C9  : num [1:3704] 21.5 20.9 20.9 22.8 NA ...
  ..$ C10 : num [1:3704] 20.7 NA NA 21.5 NA ...
  ..$ C11 : num [1:3704] 21.5 20.2 NA 22.7 20.6 ...
  ..$ C12 : num [1:3704] 20.3 NA NA 20.8 NA ...
  ..$ CS13: num [1:3704] 21.8 20.9 20.2 22.6 20.5 ...
  ..$ C14 : num [1:3704] NA NA NA 20.6 NA ...
  ..$ C15 : num [1:3704] 22 20.8 20.9 22.9 20.7 ...
  ..$ C16 : num [1:3704] NA NA NA 20.9 NA ...
 $ f_data:'data.frame': 32 obs. of  7 variables:
  ..$ Sample_ID  : chr [1:32] "H1" "H2" "H3" "H4" ...
  ..$ Sample.name: chr [1:32] "H2O" "H2M" "H4O" "H4M" ...
  ..$ Horizon    : chr [1:32] "Organic" "Mineral" "Organic" "Mineral" ...
  ..$ Treatment  : chr [1:32] "Heated" "Heated" "Heated" "Heated" ...
  ..$ T_H        : chr [1:32] "Organic_Heated" "Mineral_Heated" "Organic_Heated" "Mineral_Heated" ...
  ..$ Sample.ID  : chr [1:32] "BW.H.2.O" "BW.H.2.M" "BW.H.4.O" "BW.H.4.M" ...
  ..$ wMT        : chr [1:32] "Y" "Y" "Y" "Y" ...
 $ e_meta:'data.frame': 3704 obs. of  27 variables:
  ..$ Mass         : chr [1:3704] "200.0564425" "200.092833" "201.004077" "201.0404701" ...
  ..$ C            : int [1:3704] 8 9 7 8 7 9 8 10 9 11 ...
  ..$ H            : int [1:3704] 11 15 6 10 10 14 14 18 18 22 ...
  ..$ O            : int [1:3704] 5 4 7 6 5 5 4 4 3 3 ...
  ..$ N            : int [1:3704] 1 1 0 0 2 0 2 0 2 0 ...
  ..$ C13          : int [1:3704] 0 0 0 0 0 0 0 0 0 0 ...
  ..$ S            : int [1:3704] 0 0 0 0 0 0 0 0 0 0 ...
  ..$ P            : int [1:3704] 0 0 0 0 0 0 0 0 0 0 ...
  ..$ Error_ppm    : num [1:3704] 0.01752 -0.00742 -0.00455 -0.04212 0.02847 ...
  ..$ NeutralMass  : num [1:3704] 201 201 202 202 202 ...
  ..$ MolForm      : chr [1:3704] "C8H11O5N" "C9H15O4N" "C7H6O7" "C8H10O6" ...
  ..$ AI           : num [1:3704] 0 0 0 0 0 0 0 0 0 0 ...
  ..$ AI_Mod       : num [1:3704] 0.111 0 0.429 0.2 0 ...
  ..$ DBE          : num [1:3704] 4 3 5 4 4 3 3 2 2 1 ...
  ..$ DBE_O        : num [1:3704] -1 -1 -2 -2 -1 -2 -1 -2 -1 -2 ...
  ..$ DBE_AI       : num [1:3704] -2 -2 -2 -2 -3 -2 -3 -2 -3 -2 ...
  ..$ GFE          : num [1:3704] 53.2 73 27.7 53.2 35.9 ...
  ..$ kmass.CH2    : num [1:3704] 200 200 201 201 201 ...
  ..$ kdefect.CH2  : num [1:3704] 0.167 0.131 0.22 0.184 0.173 ...
  ..$ NOSC         : num [1:3704] 0.25 -0.444 1.143 0.25 0.857 ...
  ..$ OtoC_ratio   : num [1:3704] 0.625 0.444 1 0.75 0.714 ...
  ..$ HtoC_ratio   : num [1:3704] 1.375 1.667 0.857 1.25 1.429 ...
  ..$ NtoC_ratio   : num [1:3704] 0.125 0.111 0 0 0.286 ...
  ..$ PtoC_ratio   : num [1:3704] 0 0 0 0 0 0 0 0 0 0 ...
  ..$ NtoP_ratio   : num [1:3704] NA NA NA NA NA NA NA NA NA NA ...
  ..$ ElComposition: chr [1:3704] "CHON" "CHON" "CHO" "CHO" ...
  ..$ bs1_class    : chr [1:3704] "Lignin;Amino Sugar" "Protein" "Tannin" "Tannin" ...
 - attr(*, "cnames")=List of 29
  ..$ edata_cname     : chr "Mass"
  ..$ fdata_cname     : chr "Sample_ID"
  ..$ mass_cname      : chr "Mass"
  ..$ extraction_cname: NULL
  ..$ mf_cname        : chr "MolForm"
  ..$ c_cname         : chr "C"
  ..$ h_cname         : chr "H"
  ..$ o_cname         : chr "O"
  ..$ n_cname         : chr "N"
  ..$ s_cname         : chr "S"
  ..$ p_cname         : chr "P"
  ..$ isotopic_cname  : chr "C13"
  ..$ o2c_cname       : chr "OtoC_ratio"
  ..$ h2c_cname       : chr "HtoC_ratio"
  ..$ kmass_cname     : Named chr "kmass.CH2"
  .. ..- attr(*, "names")= chr "CH2"
  ..$ kdefect_cname   : Named chr "kdefect.CH2"
  .. ..- attr(*, "names")= chr "CH2"
  ..$ nosc_cname      : chr "NOSC"
  ..$ gfe_cname       : chr "GFE"
  ..$ mfname_cname    : NULL
  ..$ aroma_cname     : chr "AI"
  ..$ modaroma_cname  : chr "AI_Mod"
  ..$ dbe_cname       : chr "DBE"
  ..$ dbeo_cname      : chr "DBE_O"
  ..$ dbeai_cname     : chr "DBE_AI"
  ..$ elcomp_cname    : chr "ElComposition"
  ..$ n2c_cname       : chr "NtoC_ratio"
  ..$ p2c_cname       : chr "PtoC_ratio"
  ..$ n2p_cname       : chr "NtoP_ratio"
  ..$ bs1class_cname  : chr "bs1_class"
 - attr(*, "data_info")=List of 2
  ..$ data_scale     : chr "log2"
  ..$ instrument_type: chr "12T"
 - attr(*, "class")= chr [1:2] "peakData" "ftmsData"
 - attr(*, "filters")=List of 3
  ..$ massFilt    :List of 3
  .. ..$ report_text: chr "A mass filter was applied to the data, removing Masses that had a mass less than 200 or a mass greater than 900"| __truncated__
  .. ..$ threshold  : num [1:2] 200 900
  .. ..$ filtered   : chr [1:6869] "98.656748" "98.833389" "98.833638" "98.85848" ...
  ..$ moleculeFilt:List of 3
  .. ..$ report_text: chr "A molecule filter was applied to the data, removing Masses that were present in fewer than 2 samples. A total o"| __truncated__
  .. ..$ threshold  : num 2
  .. ..$ filtered   : chr [1:4258] "200.020042" "200.032592" "200.040329" "200.046097" ...
  ..$ formulaFilt :List of 3
  .. ..$ report_text: chr "A formula filter was applied to the data, removing Masses that had NoFormula assigned. A total of 612 Masses we"| __truncated__
  .. ..$ threshold  : chr "NoFormula"
  .. ..$ filtered   : chr [1:612] "200.9892554" "201.053514" "202.027165" "202.0376923" ...

peakObj

peakData object
# Peaks: 3704
# Samples: 32
Meta data columns: [Mass, C, H, O, N, C13, S, P, Error_ppm, NeutralMass, MolForm, AI, AI_Mod, DBE, DBE_O, DBE_AI, GFE, kmass.CH2, kdefect.CH2, NOSC, OtoC_ratio, HtoC_ratio, NtoC_ratio, PtoC_ratio, NtoP_ratio, ElComposition, bs1_class]
A mass filter was applied to the data, removing Masses that had a mass less than 200 or a mass greater than 900. A total of 6869 Masses were filtered out of the dataset by this filter.
A molecule filter was applied to the data, removing Masses that were present in fewer than 2 samples. A total of 4258 Masses were filtered out of the dataset by this filter.
A formula filter was applied to the data, removing Masses that had NoFormula assigned. A total of 612 Masses were filtered out of the dataset by this filter.

Binary Matrix

# Set rownames
edata <- peakObj$e_data
rownames(edata) <- as.character(edata$Mass)
edata$Mass <- NULL
# Binary matrix
edata.binary <- apply(edata,2,function(y)as.numeric(y>0))
# Na's as zero
edata.binary[is.na(edata.binary)] <- 0
# As dataframe edata.binary
edata.binary <- as.data.frame(edata.binary)
# Add Mass to binary matrix
edata.binary$Mass <- rownames(edata)
# Add metadata to binary matrix
tmp <- peakObj$e_meta[,c("Mass","bs1_class")] # save to tmp; direct merge isn't working
str(tmp)

'data.frame':   3704 obs. of  2 variables:
 $ Mass     : chr  "200.0564425" "200.092833" "201.004077" "201.0404701" ...
 $ bs1_class: chr  "Lignin;Amino Sugar" "Protein" "Tannin" "Tannin" ...

# Merge 
edata.binary <- merge(edata.binary,tmp,by="Mass",all = T)
# Remove Mass column from df
edata.binary$Mass <- NULL

# Carbon compounds #
carbon.cpds <- aggregate(.~bs1_class,data=edata.binary,FUN=sum)
# Transfor count to relative abundance #
carbon.cpds[,-1] <- decostand(carbon.cpds[,-1],method = "total",MARGIN = 2)*100
# Sum check
colSums(carbon.cpds[,-1])

  H1   H2   H3   H4   H5   H6   H7   H8   H9  H10  H11  H12  H13  H14  H15  H16   C1   C2   C3   C4   C5   C6   C7   C8   C9  C10 
 100  100  100  100  100  100  100  100  100  100  100  100  100  100  100  100  100  100  100  100  100  100  100  100  100  100 
 C11  C12 CS13  C14  C15  C16 
 100  100  100  100  100  100 

# melt df #
carbon.cpds <- melt(carbon.cpds,variable.name = "Sample_ID")

Using bs1_class as id variables

# Add sample info
carbon.cpds <- merge(carbon.cpds,peakObj$f_data,by="Sample_ID")
head(carbon.cpds)

# Color palette #
pal12 <- c("#A6CEE3", "#1F78B4", "#B2DF8A", "#33A02C", "#FB9A99","#E31A1C", "#FDBF6F", "#FF7F00", "#CAB2D6", "#6A3D9A","#FFFF99", "#B15928")
# New color palette #
pal4 <- c("#1F78B4", "#FF7F00", "#33A02C", "#E31A1C")
# Order levels
carbon.cpds$Horizon <- factor(carbon.cpds$Horizon,levels = c("Organic","Mineral"))

# Plotting
p <- ggplot(carbon.cpds,aes(x=bs1_class,y=value,fill=Treatment))+ geom_boxplot()+
      ylab("Percentage") + xlab("") + scale_fill_manual(values=pal4) + 
      theme(axis.text.x=element_text(size=12,face="bold",angle = 90,vjust = 0.7),
            axis.text.y=element_text(size=18,face="bold"),
            strip.text.x = element_text(size = 16),
            legend.text=element_text(size=18))
p + facet_grid(~Horizon)

Ordination Plot:Treatment

rm(list = ls())
load("ficusPeakObj.Robj")
peakObj <- group_designation(peakObj, main_effects=c("Treatment"))
pcoaMat <- getPrincipalCoordinates(ftmsObj = peakObj,dist_metric = "bray")

Data has already been transformed. Data is on the log2 scale. Changing to pres scale.

plotPrincipalCoordinates(pcoaMat,x=1,y=2,ftmsObj = peakObj )

minimal value for n is 3, returning requested palette with 3 different levels
minimal value for n is 3, returning requested palette with 3 different levels
minimal value for n is 3, returning requested palette with 3 different levels
minimal value for n is 3, returning requested palette with 3 different levels

Ordination Plot:Treatment x Horizon

peakObj <- group_designation(peakObj, main_effects=c("T_H"))
#getGroupDF(peakObj)
pcoaMat <- getPrincipalCoordinates(ftmsObj = peakObj,dist_metric = "bray")

Data has already been transformed. Data is on the log2 scale. Changing to pres scale.

plotPrincipalCoordinates(pcoaMat,x=1,y=2,ftmsObj = peakObj )

Venny Diagram

rm(list = ls())
# load organic
load("organicUnique.Robj")
# get values
organic.unique.df <- organic_unique$e_data
# load mineral
load("mineralUnique.Robj")
# get values
mineral.unique.df <- mineral_unique$e_data
# Organic data
All.organic <- organic.unique.df[organic.unique.df$uniqueness_nsamps %in% "Observed in Both","Mass"]
# Mineral data
All.mineral <- mineral.unique.df[mineral.unique.df$uniqueness_nsamps %in% "Observed in Both","Mass"]
# Create List
venn.data <- list(
  Control.Mineral = c(mineral.unique.df[mineral.unique.df$uniqueness_nsamps %in% "Unique to Control","Mass"],All.mineral),
  Control.Organic = c(organic.unique.df[organic.unique.df$uniqueness_nsamps %in% "Unique to Control","Mass"],All.organic),
  Heated.Mineral = c(mineral.unique.df[mineral.unique.df$uniqueness_nsamps %in% "Unique to Heated","Mass"],All.mineral),
  Heated.Organic = c(organic.unique.df[organic.unique.df$uniqueness_nsamps %in% "Unique to Heated","Mass"],All.organic)
  )

# Venn plot
ggVennDiagram(venn.data,label_alpha = 0)